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Class-Shape Transformation (CST) is a popular method for creating analytical repre- 
sentations of the surface coordinates of various components of aerospace vehicles. A wide 
variety of two- and three-dimensional shapes can be represented analytically using only 
a modest number of parameters, and the surface representation is smooth and continu- 
ous to as fine a degree as desired. This paper expands upon the original two-dimensional 
representation of airfoils to develop a generalized three-dimensional CST parametrization 
scheme that is suitable for a wider range of aircraft wings than previous formulations, in- 
cluding wings with significant non-planar shapes such as blended winglets and box wings. 
The method uses individual functions for the spanwise variation of airfoil shape, chord, 
thickness, twist, and reference axis coordinates to build up the complete wing shape. An 
alternative formulation parametrizes the slopes of the reference axis coordinates in order 
to relate the spanwise variation to the tangents of the sweep and dihedral angles. Also 
discussed are methods for fitting existing wing surface coordinates, including the use of 
piecewise equations to handle discontinuities, and mathematical formulations of geometric 
continuity constraints. A subsonic transport wing model is used as an example problem 
to illustrate the application of the methodology and to quantify the effects of piecewise 
representation and curvature constraints. 


Nomenclature 


a Bezier coefficient 

h wing span 

hj * 1 j th Bernstein basis polynomial of degree N 

c airfoil chord length 

c airfoil chord length, normalized by semispan 

C class function 

i airfoil incidence 

j index of summation 

k section index 

l index of summation 

number of terms in Bezier curve 

, N'2 class function leading-edge and trailing-edge parameters 
s spanwise arc length 

S shape function 

t airfoil maximum thickness 

x , y , z longitudinal, lateral and vertical coordinates in aircraft body axes 
x, y , z coordinates normalized by semispan 

r wing dihedral 

£ airfoil coordinate normal to chord line, normalized by chord 

£ value of £ scaled by max. thickness 

77 spanwise coordinate, normalized by semispan 

77 value of 77 normalized by local section length 
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r) value of 77 for a known data point 

A c /4 wing quarter-chord sweep 

X generic spanwise parameter 

X value of x for a known data point 

7/) airfoil chordwise coordinate, normalized by chord 

Subscripts 

L lower surface 

root value at section root 

tip value at section tip 

T trailing edge 

U upper surface 

A differential 

0 reference axis 


I. Introduction 


T he Class- Shape Transformation (CST) parametrization method ’ ’ ' has become an increasingly popu- 
lar method for creating analytical representations of the surface coordinates of various components of 
aerospace vehicles. 1,5 >< ’ CST parametrization can be used to represent a wide variety of two- and three- 
dimensional shapes with a modest number of parameters, and it provides built-in design variables for use 
in shape optimization. The use of analytical functions means that the surface representation is smooth and 
continuous to as fine a degree as desired, without the need for interpolation between discrete points. The 
analytical functions also provide the ability to directly compute surface derivatives, normals, and integrals, 
instead of approximating them using finite-difference methods. 

The basic two-dimensional CST methodology describes the normalized surface coordinates of an airfoil 
at zero incidence using analytical functions of the form C = C (VO : (0 < "0 < 1) , where if = (x — xj e ) /c is 
the normalized chordwise coordinate and C, = (z — zj e ) /c is the coordinate normal to the chord line (Fig. 1). 
The surface functions take the form of the product of a class function, (7^0 (’if), and a shape function, S (if), 
with separate functions for each of the upper and lower surfaces: 


Cu«0= c^ 2 (if)s v (if) + if( Tu 
C L (VO = c%(if)s L (if) + if{T L 


(i) 


The additional Qr terms specify the trailing-edge thickness-to-chord ratio for blunt-trailing-edge airfoils. In 
Eq. 1, the class function, , takes the form 


c$ (VO = V^ 1 (1 - V0 Wa (0 < V> < 1) (2) 

where exponents Ni and N 2 are parameters that define the type of class being represented by the function. 
Different choices for the values of these parameters can result in a wide range of basic geometry shapes 
(Table 1). The blunt-nosed airfoil class function (Ni = 0.5 and N 2 = 1) results in a symmetric airfoil with 
rounded leading edge and sharp trailing edge that is smooth and continuous everywhere in between (Fig. 2). 
In most applications, the shape function takes the form of a Bezier curve of order N: 


N ,r N 

s (VO = E aj (") r (1 - VO" = E a 3 b? (if) 


j=o 


3=0 


(0<if<l) 


(3) 


where ('J) is the binomial coefficient and the set {a-j : j = 0 . . . N} consists of N + 1 weighting coefficients 
for each of the terms in the summation. The weighting coefficients can be used as design variables to modify 
the shape of the airfoil while maintaining a smooth shape at all times. For an infinite number of terms in 
the series, the Weierstrass approximation theorem shows that there always exists a set of coefficients that 
will match any smooth shape exactly. s In addition, if the order of the series is finite, there always exists a 
sufficiently large number of terms to bound the approximation error within a desired magnitude. A Bezier 
form of shape function is assumed for all applications in this paper. 
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Figure 1. Definition of normalized airfoil coordinates 


Table 1. Sample values for coeffi- 
cients JVi and N 2 


iVi 

n 2 

Class 

0 

0 

Unit function* 

0.5 

1 

Blunt-nosed airfoil 

0.5 

0.5 

Elliptical airfoil 

1 

0 

Wedge 

0.75 

0.75 

Sears-Haack body 

1 

1 

Biconvex airfoil 


In this paper it is assumed that 0° = 1. 



Figure 2. Basic airfoil shape defined by the Kulfan blunt-nosed airfoil class function, C®' 5 
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As suggested by Kulfan, two-dimensional parametrization of the airfoil shape can be extended to three 
dimensions by representing the wing as an extrusion of parametrized airfoils along a spanwise axis. If one 
introduces an additional spanwise parameter, 77, where ?/ = 0 represents the root of the wing and rj = 1 
represents the tip, one can represent the airfoil shape anywhere along the span using an analytical function 
of the two variables if) and 77: 


Cu 


Cl (iM) 

Here, the normalized coordinate £ = -q~ c is used to allow the maximum thickness-to-chord ratio to be 
scaled independently, as described in Section III. These equations function as a transformation from the 
non-dimensional (■*/>, 77) space to the physical domain. They are controlled by a set of Nu Njj^ + Nl^Nl^ 
weighting coefficients to describe the upper- and lower-surface shapes of the airfoils for all locations on the 
wing. 

When the airfoil functions are combined with known spanwise variations of chord, twist, and leading- 
edge location, the three-dimensional wing surface can be assembled. In several studies, the planform shape 
has been specified using simple linear parametrizations for a known class of wing (e.g. a simple subsonic 
transport wing' 1 or a cranked delta wing ’). The planform shape has also been specified by defining the chord 
and leading-edge x and z coordinates using additional analytical functions of 77 = 2 y/b. 

However, these previous three-dimensional applications lack a methodology of sufficient generality to 
define all classes of wings, including those with extreme dihedral, such as blended winglets, or those for 
which the leading-edge x and z coordinates cannot be formulated as monotonic functions of y, such as box 
wings. In addition, many wings exhibit discontinuities in the spanwise variation of parameters, such as chord 
breaks, so the use of a single, continuous function to describe these parameters can be problematic. This 
paper proposes a methodology for using spanwise parametrization of the planform shape that retains sufficient 
generality to be suitable for a much wider range of extruded shapes. Section II formulates the mathematical 
representation of spanwise-varying parameters, including a method of accounting for discontinuities using 
piecewise functions in combination with geometric continuity constraints. Section III applies these concepts 
to formulate a comprehensive representation of general three-dimensional wings. Finally, Section IV applies 
the methodology to an example problem, that of a subsonic transport wing with a blended winglet. 

II. Spanwise Parametrization 

In addition to three-dimensional parametrization of the airfoils (Eq. 4) , this paper introduces a spanwise 
parametrization scheme for the variation of wing design parameters. For a given parameter, y, which 
may be any spanwise-varying property such as chord, twist, etc., one can use CST to create a functional 
representation of the property between the root and tip of a wing: 

. . N 

X (v) = C”l (77) E djbf (77) + 77X tip + (1 - v)x root (° < V < 1) (5) 

3=0 

The additional term, (1 — 77)y roo t, accounts for the fact that the property may have a non-zero value at 
the root, unlike airfoils. The following section examines the process for determining the shape function 
coefficients that best represent a known spanwise variation. 

A. Fitting CST Functions to Existing Shapes 

When using a series representation of a shape such as Eq. 5, the choice of values for the Bezier coefficients 
to produce a given shape may not be readily obvious. Instead, one would normally start with a given shape 
and solve for the values of the N + 1 coefficients aj,j = 0 ... N which best approximate it. Given M + 3 
ordered pairs of known coordinates (r)j,Xj) ,j = 0 ... M + 2, with 770 = 0 and t)m+ 2 = 1, one can insert the 
values of fjj and Xj i n turn into Eq. 5 and create a system of linear equations. Since, by identity, X r00 t = Xo 
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3 = 0 1=0 

n H n 

c nI WO H aL i.‘ b j H WO h ^ Lr> W)+ 4>Ct l 

3=0 1=0 


( o<7/><ino<7?<i ) 
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and xtip = Xm+ 2 , those terms can be moved to the right side of the equation, and the set of equations for 
each of the interior points gives M + 1 equations for N + 1 variables: 


C%l(Vi)bZ(r) 1 ) 




^ a 0 ^ 

C%(V2)bg (fh) 

Cf«2)&f(r? 2 ) 

c$(fc)$(fo) 


a\ 

. C% (ijm+i)^ (vm+i) 

(W) K (m+i) ■ 

' C%l(fjM+i)b%(fj M + 1 ) . 


\ aN / 


l ( A x)i \ 

( A X') 2 

V ( a x)m+i / 


where (Ax) j = - %Xm + 2 - (1 - Vi) Xo- 

If the number of interior points to fit is equal to the number of coefficients (M = iV), the exact solution 
for the coefficients a,j can be found through direct matrix inversion or a comparable iterative method. If 
the number of points is greater than the number of coefficients ( M > N), the coefficients can instead be 
approximated through a pseudo-inverse solution (least-squares fit). 

If the shape to be approximated is very non-linear or the interval length between points varies significantly, 
using an exact solution to Eq. 6 can result in over-fitting of the points and the introduction of spurious 
oscillations in the intervals between the known points; using a least-squares fit in these cases can help 
alleviate these problems. In general, experience has shown that it is best to use a number of known points 
that is two to three times larger than the number of coefficients to be estimated. 

Figure 3(a) gives an example of a shape that is particularly challenging: the leftmost and rightmost 
sections of the shape are purely linear, with a rapid non-linear change in \ i n the middle section (0.4 < 77 < 
0.6). Using an exact solution (M = N = 16) results in a function that does indeed pass through all the data 
points, but its behavior is significantly different than the underlying function from which the data points 
were sampled, particularly in the linear regions. Doubling the number of data points ( M = 35), for the same 
number of coefficients, greatly improves the stability of the fit. Note, however, that with M > N, the CST 
function only approximates the data points, which is especially obvious at 77 = 0.4 and 77 = 0.6. 


B. Piecewise Representation of Complex Shapes 

The use of CST functions to approximate shapes that do not deviate significantly from their class function 
usually produces acceptable results. In complex cases like the one shown previously in Fig. 3(a), however, 
a large number of coefficients and correspondingly large number of sampled data points may be required to 
reduce the amplitude of the oscillations between known points to an acceptable level. As an alternative, the 
shape may be represented in a piecewise manner by defining functions on two or more adjoining sub-intervals 
and using an appropriate number of Bezier coefficients for each. If the curve is divided at P + 2 break points 
Vo,Vi,V2,---Vp+i, where 770 = 0 and 77 P +1 = 1, then the equation for the generic parameter x becomes a 
combination of the equations over the P + 1 intervals: 



( N 0 

CZ (%) E a j 0 6 f° (770) 

9=0 

+ 

? 7 oXtip 0 

+ 

(1 Vo) Xrooto 

(0 < 77 < 771) 


x 0?) = < 

Ni 

Oh) E Oh) 

3=0 

+ 

ViXti Pl 

+ 

(1 — Vi) XrooU 

(jli < V < V2) 

(7) 


Np 

c£ (Vp) E a JP b?r (77 p) 
l 3= 0 

+ 

t 1 pX ti Pp 

+ 

(1 — Vp) Xrootp 

{jlP < V < 1 ) 



where 77 ,. = n nk . 

Vk+l-Vk 

Using this piecewise definition of the spanwise function, Eq. 6 may be modified to create a block diagonal 
matrix equation. In this case there are P + 1 sets of ordered pairs, with the fcth set (0 < k < P) having 
M k + 3 ordered pairs of coordinates (% fc , XjJ, j = 0 . . . M k + 2 , with T) 0k = r] k and V Mk + 2 k = Vk+ 1 - Again, by 
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(a) Single-section fit with two data set sizes 



(b) Single-section vs. piecewise fit 


Figure 3. Sample CST curve fits using a C\ class function and Bezier-polynomial shape function(s) with 16 
coefficients total 
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identity Xroot k 
variables: 


where 


p p 

\o k and Xtip = xm+ 2 k , so the system consists of (^k + 1) equations for Y2 (-Nfc + 1) 
Vk k = 0 k = 0 


\Ck\ = 


1 

0 

0 


( [«o] \ 


f bco] \ 

0 [Cl] ••• 0 


M 

= 

[xi] 

0 0 ■■■ [Cp\ _ 


V [ap] ) 


l \xp] / 


c nI (vi k )b% k (fi lk ) ••• c$(f) u )b%;(fi lh ) 

-^nI ( ? ?M fc + lfc )b jiM k +l k ) ' ' ' (^Mfe+lfc) (^M fc + l fc ) -l 

( a 0 k \ 


[afc] = 


ai 


k 


( 8 ) 


bck] 


\ a N k ) 

( ( A x)i fc \ 

( A *k 


V ( A X) M k + l k / 


while k = JI+,-4 and ( A % = Xj k ~ f) jk XM+ 2 k - (1 - fj jh ) X0 k - 

Figure 3(b) illustrates how fitting of the previous example can be greatly improved through the use of 
piecewise functions. For the piecewise case, the number of interior data points available for curve fitting 
is reduced by two because the data points at rj = 0.4 and rj = 0.6 now become end points for multiple 
segments; additionally, the number of Bezier coefficients is reduced by two because of the need to specify 
the function values at the joints (x^p = Xrooti = d and Xtip = Xroot 2 = !)• Nevertheless, the piecewise 
parametrization scheme exhibits improved performance for the same number of total data points and fitting 
coefficients. All of the Bezier coefficients are used for the middle segment, while none are required for the 
left and right segments since those regions are linear. 


C. Geometric Continuity Constraints 

The single-section parametrization scheme of Eq. 5 is continuous, as are all of its derivatives. However, 
a piecewise representation generally is not continuous across the section breaks, and its use may introduce 
discontinuities even when the sampled data points come from a continuous curve. There are cases in which it 
is desirable for the shape to exhibit a certain degree of continuity in order to produce acceptable performance, 
such as low drag, resistance to boundary-layer separation, etc. If two curves are joined, they are said to 
have zero-order geometric continuity , G°, at the joint if their end points exactly coincide. In addition, the 
curves have first-order geometric continuity , G 1 , if their tangents point in the same direction at the joint; 
and they have second-order geometric continuity , G 2 , if their centers of curvature at the joint coincide. In 
the following sections, the conditions for enforcing geometric zero-, first- and second-order continuity across 
section breaks of a piecewise parametrization scheme are derived. The continuity conditions are substituted 
into Eq. 8 to derive a modified form of the linear set of equations. 

The behavior of the class function, G^ 1 ( 77 ), at its end points (rj = 0 and g = 1) can vary greatly depending 
on the values of Ni and _/V 2 ; the behavior at the end points in turn affects the form of the continuity constraint 
equations. For the biconvex class function, C\ ( 77 ), the values of the function and its derivatives at the end 
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points are 


are 


Cl (0) = 

0 , 

C\{ 1) 

= 0 


C'\ (0) = 

1 . 

C'\{ 1) 

= -1 

(9) 

C"\ (0) = 

- 2 , 

C"\{ 1) 

= -2 


b Co ( 77 ), the 

values of the 

function and its derivatives at the end points 

C 0 °( 0) = 

1, 

Cg(l) 

= 1 


II 

O O 

b 

0 , 

C'o (1) 

= 0 

(10) 

II 

0 0 

b 

0 , 

C"o(l) 

= 0 



The following discussion is limited to the derivation of continuity constraints for these two forms of class 
function, in addition to the blunt-nosed airfoil class function, C J’ 5 ( 77 ). 


1. Zero- Order Continuity 

For the piecewise parametrization of generic parameter \ (Eq. 7), the values of the function at the root and 
tip of section k are 

Xk {ilk) = C N ^{0)ao k + Xrootie 
Xk{Vk+ 1 ) = C%l{l)a Nk + y tipfc 

If the tip of segment k is joined to the root of segment k + 1, enforcing G° continuity at the joint requires 
that Xk (jlk+i ) = Xk + 1 (jlk+i )• For the biconvex class function, the value is zero at the end points so ensuring 
G° continuity requires only that 

Xroot fc+ i = Xtip fe (12) 

For a pure Bezier parametrization, there is an additional requirement that 


a 0 fc+ i — ®N k 


(13) 


2. First- Order Continuity 

The equations for the derivative of the CST parametrization of section k at the root and tip are 

Xk {Vk) = (0) a 0fc + NkCff* (0) (— a 0k + ai fc )+ A\ k 

x'k iVk+l) = C nI (1) a N k + NkC^l (1) ( — 0-N-\ k + «ATfc)+ AXk 

where A Xk = X'tip — X'root fc - For C 1 continuity at the joint between sections k and k + 1, the slopes of the 
two curves must be equal at that point: 


Xk (Vk+i) = Xfc+i ( ? ?fc+i) 
Arj k Arik+i 


(15) 


where Ai] k = {Vk + 1 — 77 fc). For the biconvex class function, one can combine the G° continuity condition 
(Eq. 12) with Eqns. 9, 14, and 15, then solve for ao k+1 as a function of a,N k : 


a o k+1 


Aiik+i 
A rj k 


a N k + 


Aiik+i 

Arjk 


Axk - Axk+i 


(16) 


For a pure Bezier parametrization, one instead uses Eqns. 10, 13, 14, and 15 to solve for ai k+1 as a function 
of a Nk and a N _ lk \ 


a i fc+ i - 


Arjk+i N k ( Ar^k+i N k \ 1 

+ l 1 + + JVwT 


(^f A “ - A “ +1 ) 


(17) 


For a blunt-nosed airfoil {Ni = 0.5), the derivatives of the upper- and lower-surface functions are positive 
and negative infinity, respectively, at the leading edge, so first-order continuity is automatically maintained 
between the upper and lower surfaces. 
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3. Second- Order Continuity 

The equations for the second derivative of the CST parametrization of section k at the root and tip are 


x! k (Vk) =C" N n \ (0) a 0k + 2 N k C'%l (0) (-a 0fe +a lh ) + 
N k (N k -l)C$(0) (a 0k -2a lk +a 2k ) 

Xk (Vk+i) =C" N \ (1) a,N k + 2 NhC' N l (1) (— ajv- i k + a,N k ) + 

Nk (. Nk — 1) Cff* (1) ( a N-2 k — 2aAr-i fc + aw k ) 


For G 2 continuity at the joint between sections k and k + 1, the centers of curvature for the two curves must 
be the same. As a necessary condition, the radii of curvature of the two curves must be equal, and since it 
is already established from G 1 continuity that the magnitudes of the first derivatives are equal, this means 
that the magnitude of the second derivatives of the two curves must also be equal at the joint: 


Xk (%+ 1) 


Xfc+i (Vk+ 1) 

(A rjk) 2 


(A?7fc +1 ) 2 


(19) 


For section k connected at the tip to section k + 1 at the root, the sign of the second derivatives must be the 
same. For two sections connected at their roots (such as two airfoil surfaces connected at the leading edge), 
the second derivatives would be of equal magnitude but opposite sign. 

For the biconvex class function, one can combine the G 1 continuity condition (Eq. 16) with Eqns. 9, 18, 
and 19, then solve for a± k+1 : 


°ifc+i — 


f Ap k+1 \ * N k 


\ A rjk ) N k+ 1 
Afc+i + 1 / Ar ) k+ 1 


a N- i fc — 


( ^Vk+i \ Af k + 1 A-q k+ i N k+1 + 1 


V 

A Xk - Ax fc+ i 


N k +i A?y k N k+ 1 


a N k 


N k+ 1 V Ik 

For a pure Bezier parametrization, one instead uses Eqns. 10, 13, 17, 18, and 19 to solve for a 2k+1 ' 

Ari k+ i N k - 1 


f Ar] k+1 \ N k N k -1 0 A?fc+i N k 

a 2*+i - i 7 a ^-2 fc ^ 2 - 


\ A?7fc J N k+ 1 N k+ 1 — 1 


1 + 


Ai-f k N k+ 1 [ A?7fe N k+ 1 - 1 


a N-l k + 


f A?7fc +] 

\ A rj k J N k+ i N k+ 1 — 1 


N k Af k — 1 Ar] k+ i N k 


A N k+ 1 


1 


a Nk 


f A r/ k+ 


N, 


k + 1 


V A77 k 


-A Xk - Axfc+ij 


(20) 


( 21 ) 


For the the blunt-nosed airfoil class, the second derivative of both the upper- and lower-surface curves 
is infinite at the leading edge, but Kulfan has shown that the radius of curvature is finite and that G 2 
continuity can be enforced simply by setting ao fc+1 = — ao k . 


III. CST Parametrization of an Aircraft Wing 

This section derives a complete set of functions to represent a wing using CST parametrization in the 
form of Eq. 5. One can define a reference axis , which is the locus of the quarter-chord points of the airfoils 
from root (?? = 0) to tip ( r] = 1). The coordinates along this axis are (xq, 2/0 i 2 o)> an d the spanwise distance, 
s, is the arc length along a projection of the reference axis onto the x = 0 plane: 

s (v) = \]vl + (22) 

The reference span, b, is twice the spanwise distance at the tip ( b = 2s \ v= i ). 

If one considers the wing to be a continuous extrusion of airfoils, one can define parametrizations that 
define the variation of chord, c, tliickness-to-chord ratio, t/c, and incidence, i, of the airfoil along the reference 
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r ? c tip+ 


(1 ??)c ro 0 t 


N c 


-(v) = 


- (v) = 

C 


i (v) = 


c nI fa) & f c fa)+ 

J= 0 

fa)XX' 6 f‘ (^)+ »?(£) + 

\ c /tip \ c /root 

Ni 

C nI fa) Y1 a h b T fa)+ 7 ?*tip+ fa - 7 ?)®root 

i=o 


(0 < ?? < 1) 


(23) 


where c = -P-. These quantities are familiar to the aircraft designer and their spanwise variations follow 
directly from non-dimensional wing design parameters such as aspect ratio, taper ratio, washout, etc. They 
represent scaling and rotation operations applied to the normalized airfoil surface obtained from Eqs. 4. Any 
or all of the sets of Bezier coefficients ({a Cj }, {«q } or {a^}) ma Y be used as design variables in shape 
optimization studies. Using separate equations for these parameters means that an optimization study could 
be performed on one of the quantities of interest independently; for example, the coefficients of the incidence 
equation could be chosen to match a desired spanwise lift distribution, all while keeping the airfoil, chord, 
and thickness distributions constant. 

One or more of these parametrizations could also be formulated as piecewise equations in the form of 
Eq. 7, with or without continuity constraints. The form of each equation can be specified separately, so that 
the chord equation might be defined in a piecewise manner with first-order discontinuity at a chord break; 
whereas the thickness and incidence equations might, at the same time, be defined as continuous across the 
break, or even as single-section parametrizations across the entire span. The use of separate equations give 
the designer a great deal of flexibility to choose the best form of equation for each parameter. 


A. Parametrization of Reference Axis Coordinates 

In addition to the distribution of the chord, thickness, and incidence, one also needs to know the physical 
locations of the airfoils along the span before one can assemble the full wing shape. One method of specifying 
these reference coordinates is to parametrize them directly using separate functions of rj , as follows: 


Xo (r/) = 

c nI b/)X/'-'vb " M + 

j=o 

^ti P + 

(1 V) x TOOt 



Vo (v) = 

Ny 

c nI fa)+ 

i=o 

wtip+ 

(! - v)Vroot 

(0 < r) < 1) 

(24) 

(v) = 

N z 

c$fa)£°*i 6 f'fa)+ 

^tip+ 

(1 — ??)~root 




where Xq = 2 ) : ° , y 0 = and Zq = ~l" ■ This direct parametrization scheme is simple and straightforward; 
on the other hand, the functions themselves may not offer insight to the aircraft designer, who is used to 
working with such wing shape parameters as sweep and dihedral. 

It is apparent that the quarter-chord sweep, A|, and dihedral, T, are related to the derivatives of the 
reference axis coordinates as follows: 


=tan 1 

4 

T =tan _1 



= tan 1 


/ d2o/ds\ 

Vdyo/ds/ 


(25) 


Therefore, it may be preferable to use an alternative 


form of the reference axis equations, written in terms 
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of the derivatives: 


(^ 7 ) M = c nI (v) £ a *> b f fa) 

v 7 i=o 

N 

(^f) (??) =c ^ (??) i a % 6 f a w 

(77) (??) =c $ w 

V 7 J=0 


/ dx 0 \ 

V ds / tip 


+ (1-7?) 


( dx( A 

V ds / roa t; 


+?7 


( d^o \ 

V ds / tip 


+(!-??) 


( d^o \ 

V ds / root 


/ dz 0 \ 

V ds / tip 


+ ( 1 - 7 ?) 


dzo\ 
ds / root 


(0 < 7? < 1) 


(26) 


In this alternative form, the reference axis is derived by integrating the derivative equations in the spanwise 
direction: 


x o ( t ?) 



d7? 



20 (t?) = 




d7? 


(27) 


Marshall has shown that CST parametrizations of all class functions with integer class parameters Ni and 
N 2 , as well as the blunt-nosed airfoil class, can be expressed as pure Bezier curves; and since Bezier curves 
have closed-form integrals, 8 this means that the integrals for these classes retain the analytical nature of the 
CST parametrization scheme. 

Although they are only directly related to the tangents of the sweep and dihedral, the alternative formu- 
lation of Eq. 26 nonetheless can be a more intuitive form for the aircraft designer. In addition, the values 
for the sweep and dihedral are often constant across the span so these equations can be simple in their form. 
Depending on the needs of the designer, either Eqs. 24 or Eqs. 26 may be used to define the coordinates of 
the reference axis. 


B. Assembly of the Three-Dimensional Wing Surface 

Using the full suite of equations (Eqs. 4 and 23), plus either Eqs. 24 or 26, the three-dimensional wing surface 
is assembled as follows: 


1. For any given point in untransformed space, (ip, rf), the normalized coordinates of the airfoil upper and 
lower surfaces are determined from Qu (ip,rj) and C, L (ip, ??). 


2. The chordwise coordinate of the airfoil is scaled by |c (rf) and the normal coordinate is scaled by 

|c ( t?) I (??). 


3. The normalized coordinates of the spanwise reference point (xq, y 0 , zq) are determined and scaled by 
|, and the dihedral of the section is determined from the derivatives of the spanwise equations using 
the second of Eqs. 25. 


4. The scaled section is rotated about its quarter-chord by the incidence angle and then about the x-axis 
by the dihedral angle, and finally it is translated so that its quarter-chord point coincides with the 
reference point. 
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In matrix notation, the above transformation steps can be expressed as 


/ 

* \ b ( xo \ b 

’ 1 

0 

0 


y\= -J + 5* 

0 

cosT 

sinT 

V 

z) \zo) 

0 

— sinT 

cosT 


cos i 0 — sin i 

0 10 
sin i 0 cos i 


0 

(£)c 


( Xo ) 


b_ 

cos i 

0 

— sini 


( 


y 0 

+ 

2° 

sin r sin i 

cosT 

sin T cos i 



° 

U ) 



cos F sin i 

— sinT 

cos T cos i 


\ 

(!)c / 


IV. Subsonic Transport Wing Example 


(28) 


As an example case, consider a notional subsonic transport wing with a blended winglet (Fig. 4). The 
wing uses a constant 10% thick supercritical airfoil and has no twist. The quarter-chord sweep and dihedral 
are constant from the centerline to the base of the winglet, and the chord exhibits a linear variation from the 
centerline to approximately 28% of the semispan, and linear variation from there to the base of the winglet. 
The combination of linear variations in the chord, sweep, and dihedral out to the tip of the main wing, 
followed by a dramatic change in the shape of the wing in the winglet, poses a challenge when attempting 
to fit this wing using spanwise CST functions. 



Figure 4. OpenVSP model of a transport wing 

The transport wing was modeled in OpenVSP 2 and exported to a set of discrete points with a grid 
resolution of 193 cross sections and 101 points per cross section. To test the ability of the different CST 
parametrization schemes to approximate the surface shape of the transport wing, three types of parametriza- 
tion were tested for this example: single-section, piecewise, and piecewise with continuity constraints. 

A. Single-Section Parametrization 

The first example case attempts to fit the transport wing geometry with a single-section parametrization. 
Sample CST curve fits for chord and dzo/ds are shown in Fig. 5(a). The distribution for chord is linear for 
the whole wing, including the winglet, whereas the distribution of dzo/ds is only linear from the wing root to 
the tip of the main wing (root of the winglet), at which point it increases dramatically as the winglet curves 
upward. Consequently, the single-section curve for the chord distribution matches the actual distribution 
well. On the other hand, a single-section CST curve for d^o/ds has difficulty fitting both the constant main- 
wing section and the greatly-varying winglet section, resulting in a compromise between the two regions and 
a number of extraneous oscillations in the winglet region. 

Figure 6(a) shows contours of fitting error superimposed over the original geometry. At each point on 
the geometry, the fitting error represents the distance between the original point and the estimated point 
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Figure 5. CST curves for chord and dzo/ds compared to original transport wing 
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(a) Single-section fit 
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(b) Piecewise fit 
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(c) Piecewise fit with continuity constraints 


Figure 6. Fitting error (relative to semispan) for single-section CST parametrization of transport wing 


14 of 16 


American Institute of Aeronautics and Astronautics 


from the CST parametrization at the same fractional spanwise and cliordwise location; the error value is 
expressed as a fraction of the wing semi-span. A skeletal representation of the CST wing is also shown on 
Fig. 6(a) to illustrate visually where it differs from the original wing surface. The difficulty in matching the 
dzo/ds distribution (and also the dyo/ds distribution, not shown) results in fitting errors as high as 1.5% of 
the semispan along the winglet surface. 

B. Piecewise Parametrization 

The second example uses piecewise parametrization to improve the approximation of the wing surface. This 
parametrization scheme uses separate curves for the inboard and outboard sections of the main wing on 
either side of the chord break, plus a third section for the winglet. For this case, only G° continuity is 
maintained across each of the section boundaries. The plot of spanwise properties (Fig. 5(b)) shows that the 
use of piecewise equations solves the problem of oscillations in the curve fit seen previously. 

Figure 6(b) shows the contours of fitting error plotted on the original surface and a skeletal outline of the 
surface of the CST wing. Note that the magnitude of the contour scale is smaller by an order of magnitude 
than the scale in Fig. 6(a), showing that the maximum fitting error is less than one tenth of the maximum 
error seen in the single-section parametrization. 

C. Piecewise Parametrization with Continuity Constraints 

The final example adds continuity constraints to the piecewise parametrization of the previous example. The 
additional constraints are G 2 continuity between the upper and lower surfaces of the airfoil at the nose, plus 
G 1 continuity of the reference axis slopes (d^o/ds, dyo/ds, and dz 0 /ds) at the root of the winglet, which is 
in effect a second-order continuity condition because the constraint is imposed on the slope of the reference 
axis. 

As can be seen in Fig. 5(c), the addition of the G 1 continuity constraint to the CST function for dzo/ds 
makes an imperceptible change to the shape of the function in the vicinity of the winglet root, relative to 
the unconstrained parametrization. This is because the dihedral changes rapidly in this region and a small 
initial error is quickly corrected. 

Figure 6(c) shows that the additional continuity constraints result in a very small increase in the fitting 
error at the winglet root. This is due to the necessary compromise between minimizing the approximation 
error to the actual surface points, and matching the desired slope at the joint exactly. The original OpenVSP 
model has only zero-order continuity of the dihedral at the winglet root, so imposing first-order continuity 
slightly degrades the fidelity to the original surface coordinates. 

V. Conclusion 

This paper has laid out a comprehensive method for using three-dimensional CST parametrization to 
represent aircraft wings, and other lifting surfaces. The parameters chosen are familiar to the aircraft designer 
for describing wings, and the independent parametrization means that parameters with linear variation (such 
as chord in the transport wing example) can be represented with only a few values, whereas parameters with 
highly non-linear variation (such as dihedral in the example) can be represented with a larger number of 
values. By representing these parameters as a function of spanwise distance (s), rather than of transverse 
distance (y) as in previous studies, the three-dimensional formulation laid out in this paper is suitable for a 
wider variety of wing shapes than previous formulations — particularly for wings with significant non-planar 
components such as blended winglets and box wings. With minor modifications, the methodology could also 
be used to describe other shapes, such as fuselages and nacelles. 

Single-section parametrization can be effective in approximating the surface of a simple trapezoidal 
wing with modest variations in chord, sweep and dihedral. For more challenging cases, such as the blended 
winglet example shown here, or for advanced configurations such as a blended wing body, the use of piecewise 
parametrization can be crucial for accurately describing the surface. 

However, the use of piecewise parametrization can introduce discontinuities into what were originally fully 
smooth and continuous functions, which is a concern that has not been addressed in previous studies. This 
paper introduces modifications that can be made to the piecewise representation for enforcing continuity 
between segments. On the one hand, the use of continuity constraints tends to slightly increase the fitting 
error near the section breaks relative to the unconstrained case; on the other hand, enforcing continuity can 
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maintain certain characteristics of the wing that can be essential to its performance. For example, first- 
or even second-order spanwise continuity of the surface could help maintain laminar flow under conditions 
of modest spanwise flow. Additionally, second-order continuity at the airfoil leading edge can help reduce 
drag in off-design conditions. Enforcing continuity constraints can be especially important when the CST 
parameters are used to modify the wing surface during an optimization process, to stop what was originally 
a mild discontinuity from becoming a severe one as the optimization proceeds. 
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